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Abstract. 

We show how various mathematical formalisms, specifically the catastrophe formalism 
and group theory, aid in the study of relevant systems in quantum optics. 

We describe the phase transition of the Dicke model for a finite number N of atoms, 
via 3 different methods, which lead to universal parametric curves for the expectation 
value of the first quadrature of the electromagnetic field and the expectation value of 
the number operator, as functions of the atomic relative population. These are valid 
for all values of the matter-field coupling parameter, and valid for both the ground and 
first-excited states. 

Using these mathematical tools, the critical value of the atom-field coupling 
parameter is found as a function of the number of atoms, from which its critical 
exponent is derived. 



PACS numbers: 42.50.Ct, 03.65.Fd, 64.70.Tg 
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1. Introduction 

With the ability to manipulate single atoms and photons in a cavity came a renewed 
interest in the models that describe their behaviour. The Hamiltonian for matter and 
field interactions can be written in two equivalent forms pQ, viz., the one associated to 
the multipolar expansion 

= E (£ + "(r.)) - e ? (E . ?.) + l I (B> + B 2 ) *f (1) 

where the first term describes the atomic field, the second the dipole interaction with 
the electric field, and the third the radiation field itself; and the one related to the 
radiation gauge 

H R =f- + V(r a ) -J-Y, (A ■ A ) + ^—2$ + r I + B 2 ) (2) 
2m mCa 2mc z on J 

where the long wavelength approximation is considered, and which contains a so-called 
diamagnetic term quadratic in the electromagnetic vector potential A. They are both 
related via the (unitary) gauge transformation (TJ |2] 

e N 

U = exp[i r y £rs-A], (3) 

nC s=l 

i.e., both Hamiltonians give rise to the same physical predictions when treated correctly. 

When approximating to 2-level systems, the multipole Hamiltonian gives rise to 
the well known Dicke model (DM) [3]. This describes the interaction of a single mode 
quantized radiation field with a sample of N two-level atoms, located inside an optical 
cavity, in the dipolar approximation (i.e., located within a distance smaller than the 
wavelength of the radiation), and has the form 

1 j. 1 T / + T / + 

H D = —hu)p(va + —uaJ z -\ — t= — ( a J- + aJ+)-\ — ;= — (a 1 J + + aJ) (4) 
N N JNN \/NN 



where uf is the field frequency, TiCja is the atomic energy separation, a) and a denote 
the one-mode creation and annihilation photon operators respectively, J z the atomic 
relative population operator, and J± the atomic transition operators. We have taken 
H — 1. The interaction Hamiltonian (third and fourth terms on the right hand side) 
contains, respectively, the rotating and counter-rotating terms, with coupling constant 7. 
Note that we have divided the Hamiltonian by the total number of particles N, having in 
this way an intensive Hamiltonian operator. This is convenient for the thermodynamic 
limit is found by just taking iV — > 00, though our interest is to study the behaviour at 
finite N. 

The approximation to 2-level systems of the radiation Hamiltonian Hr, however, 
yields a different Hamiltonian Hr21 which contains the diamagnetic term and which 
is not gauge-equivalent to Hp. This is summarised in the following non- commutative 
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diagram 
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where the downward dashed arrow signifies that no gauge transformation exists between 
the two. 

We stress here the distinction between H D and H R2 £ because there has been 
confusion in the literature. An important feature of the Dicke Hamiltonian Hd is 
the presence of a phase transition from the normal to a super-radiant behaviour jj]. 
This is a collective effect involving all N atoms in the sample, where the decay 
rate is proportional to N 2 instead of N, the expected result for independent atom 
emission. This phenomenon was observed experimentally in optically pumped HF 
gas [5]. However, the transition from the normal to the super-radiant regime has been 
much debated in the literature, mainly because when using H R2 e the Thomas- Reiche- 
Kuhn sum rule would place contradictory bounds to the parameters of the model (cf. 
e.g. [6]), and because it was believed that gauge invariance requires the presence of the 
diamagnetic term [7]. (The essence behind these arguments lies in the fact that the 
coupling strength 7 is much smaller than the atomic level separation ua, for optical 
systems.) At variance with this common belief, it has more recently been argued 
that the Dicke model does give a consistent description of the interaction of a one- 
mode light field with the internal excitation of atoms inside a cavity [8]. Furthermore, 
recent experimental results indicate that it can actually be observed (9J [TO] by using 
Raman transitions between ground states in an atomic ensemble, thus circumventing 
the obstacle 7 <C uja- 

In general, if we write Hr> = Hq + V, with Uq = exp[— (i/h) Ho t] the 
evolution operator for H , then, since UqClUo = a exp[— i{up/N) t] and UqJ + U = 
J + exp[—i(u}A/N)t], we have 



V I = U j VU 



I a T e n J_ e n + a e * J + e * 



+ 



j\[i/2 \ 

= J^m ( atj - e ^ ( " F "" A) * + a J+e-M"*-*^) + 

When uf and Co a are close to one another (as when in or close to resonance), and 
of the order of optical frequencies (~ 10 15 s _1 ), the terms of the form e ±l ( w F+o; A )t var y 
rapidly in a time-scale of 1/up with an average close to zero, and may be neglected. 
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This is the well-studied Tavis-Cummings [32] model (TCM) (cf. e.g. [El El [15]), and 
the approximation is called the rotating wave approximation. 

Figure [I] shows this very clearly [TT] , where the time evolution of the excited state 
|e) of one atom has been plotted for Cja = = £l (left), and uja = = 19 ft (right), 
with £1 = 7/(iV 3 / 2 fi) the Rabi frequency of the atom; here, the time scale has been taken 
in units of 1/Q. As Cja increases, the evolution approaches the harmonic oscillations 
that can be obtained analytically in this approximation. 

tr(p|e)<e|) tr(p|e>(e|) 




Figure 1. Evolution of an excited atom. The plots show the probability of finding 
the atom in its excited state, as a function of time, without using the rotating wave 
approximation. When Co a is small (left, Cja = O) high frequency terms deviate 
the result from the harmonic behaviour obtained analytically in the rotating wave 
approximation. These terms tend to disappear when Cja is large (right, Cja = 19 0). 

The TCM model has more recently been generalised [16], and considered to describe 
cavity QED with Bose-Einstein condensates [TYJ. The simplest case iV = 1 in the 
rotating-wave approximation is known as the Jaynes-Cummings model (JCM) [18]. 
Many theoretical predictions of this latter model, such as the existence of collapse and 
revivals in the Rabi oscillations [12], the formation of macroscopic quantum states, or 
measures of entanglement associated with spin-squeezed states, have been confirmed, 
and many experimental studies of Rydberg atoms with very large principal quantum 
number within single-mode cavities have been realised [20]. Together, the JCM and the 
TCM have served as a guide to understand several quantum optics phenomena. 

In this work we will consider the full Dicke Model, with rotating and counter- 
rotating terms, and for simplicity we define to a = — , 1 = ^ = T- We shall also 
measure frequencies in units of the field frequency, i.e. ujp = 1, so that 

H = I fl t fl + -j^Jz + + a){J + + J_) (6) 

and ua, 7, J+, J-, J z , Jx, J y are all adimensional. The coupling parameter 7 is given, 
in terms of physical quantities of the system, as 7 = yj (27ig/hup) d^a • ep, where g 
denotes the density of atoms in the volume of quantisation, d\, a is the excitation matrix 
element of the electric dipole operator of a single atom, and e*p the polarisation vector. 

We will show how to calculate the phase transition of the Dicke model for a 
finite number iV of atoms, via 3 different methods: i) numerical diagonalisation of 
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the Hamiltonian and the use of the fidelity between neighbouring states; ii) variational 
test states that are a direct product of coherent Heisenberg-Weyl (HW (l))-st&tes (for 
the electromagnetic field), and S f i7(2)-states (for the atomic field) [SUES]; and iii) the 
use of projection operators on the coherent states in (ii) to obtain states which obey 
the parity symmetry in the total excitation number, present in the Hamiltonian. It will 
be seen that the symmetry- adapted states (iii) constitute a much better approximation 
(than the simple coherent states) to the exact quantum ground and first-excited states, 
apart from offering the localisation of a precursor of the quantum phase transition for 
finite N (usually referred to in the literature as a quantum phase transition in a finite 
system [25]). 

It is rather surprising that the three different approaches mentioned above lead to 
universal parametric curves for the expectation values of both the first quadrature of the 
electromagnetic field, and the atomic relative population, as functions of the atom-field 
coupling parameter. These are valid for the ground- and first-excited states, and are 
also presented. 

Finally, using these procedures, the critical value of the atom-field coupling 
parameter j c may be found as a function of the number of atoms N, from which its 
critical exponent may be derived. 



2. Coherent States 

The systems described by Eq.([6]) are non-integrable and must be solved via numerical 
diagonalisation for specific scenarios, but good semiclassical approximations that mimic 
well the exact quantum behaviour have been constructed [211 E2 12"%] . 

One of these is to take as a test state a direct product of coherent ifW / (l)-states 
\a) (for the electromagnetic field), and 5'L r (2)-states \Q (for the atomic field): 

-H 2 / 2 oo +j ( u / oi \ 1/2 ] 
(l + |C| 2 ) ^om=-, [V^.\3+mJ J 

here v represents the eigenvalues of the photon number operator a^a, m those of the 
atomic relative population operator J z , and j(j + 1) those of the squared total angular 
momentum operator J 2 . We can calculate the expectation values of the field and matter 
observables by using the representation of the angular momentum and Weyl generators 
with respect to the tensorial product \a) <8> \() [2S]. Acting with the operator J z on an 
arbitrary state and projecting it onto the SU(2) coherent state \(), we obtain 

where we have used e^* J ~ = J_e''* J -, and written \(} = (1 + \C\ 2 ) J \C) f° r t ne 



Mathematical Methods in Quantum Optics: the Dicke Model 



6 



unnormalised state. The second term in parenthesis may be calculated as the operator 
1 r 9 jCC d 

Fit* n _L IAI2V7 + 1 ^ 



:i + ici 2 )^ 



(1 + |C| 2 )^ * d(* (1 + Kl 2 )^ 1 * d(* 
so that, substituting above, we get 

<ci^w = (-r^ + c-^)«w 

i.e., the representation for J z in the basis of coherent states is 

Similarly, for the other operators we have 

j_C d_ 
1 + |C| 2 + d(* 

d a 

a -> t; 1 — 

da* 2 

a) ^a* (10) 

The energy surface is then easily found by taking the expectation value of the Dicke 
Hamiltonian with respect to \a) £g> \(), and the best variational approximation to the 
ground state energy of the system, as well as its corresponding eigenstate, are obtained 
from it. In our case the energy surface has the form 

H(a, C) = H ® (CI H \a) g> |0 

1 f 2 . 2\ U A „ , 7 



11T (f> + <? — cos6M — i=q sm.9 cos0 (11) 

2iv v ' 2 

where we have defined a = -^(q + ip), C — tan f|) exp(i0), and where (g, p) 
correspond to the expectation values of the radiation field quadratures, and (9, 0) 
determines a point on the Bloch sphere. 

Critical points which minimise the energy surface are obtained via the catastrophe 
formalism [26] 

9 C = 0, q c = 0, p c = for | | < 7 C 

9 C = arccos(7 c /7) 2 , g c = icichY cos C , p c = 0, C = 0, 7r 

for |7| > 7 C 



where 7 C = ^oj a /A. This defines the separatrix of the phase transition, as shown in 
Fig. [2] Note that, for the CS used, it is independent of N(= 2j). 

Figure [3] (left) shows the critical value for the (normalised) field quadrature q c /y/N 
as a function of the interaction strength 7. The presence of a phase transition from the 



Mathematical Methods in Quantum Optics: the Dicke Model 



7 







15 








normal 


10 


N 






region 


5 




super- nidi ant 
region 


-2 


-1 


-5 
-10 
-15 


l 


2 



Figure 2. Separatrix for the normal and super-radiant regions given by the critical 
points obtained from the coherent states. 




Figure 3. Left: Critical value of as a function of the interaction strength 7. Right: 
q c /\/N versus 9 C for all 7 and any number of atoms. In both curves oj a = 1, 4>c = 0. 



normal (where this quantity vanishes) to the super-radiant phase is evident. From the 
expressions for the critical points we obtain, in the super-radiant region, 

q c , — sin# c 

-j= = -^J A y= =g= cos <j) c (12) 

VA^ v 2 cos# c 

Figure [3] (right) shows this relationship; the normal regime is described only by the 
origin (q c = 0, 9 C — 0). The rest of the curve belongs to the super-radiant regime. Note 
that this curve is valid for all \ j\ > 7 C (as this parameter drops out from the equation 
above). It is also valid for any number N of atoms. We shall see in what follows that 
this curve is universal. 

A natural variable to define is x — 7/70 and the energy at the minima take the 
values 

-^normal — — 2A^7 C , -^super-radiant — ~^1c x y- + x ) (13) 
One may study the statistical properties of the variational states by calculating the 
expectation values of linear matter and field observables with respect to the CS, at the 
minima, as well as their fluctuations [22]. For instance, for the eigenvalues A of the 
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excitation number operator A = a) a + J z + yJ 2 + l/4 — 1/2, we have 

N 



■^normal ^ ' "^super-radiant 2 



x" 2 + 2 7c 2 x 2 (l -x" 4 )] (14) 



The fluctuations AA are zero in the normal region, and in the super-radiant phase take 
the form 



AA = y|Q+2 7c 2 x2)(l-x-^) (15) 

A complete table for the expectation values and fluctuations of the most important 
matter and field observables is given in [22]. Many of these resemble very well the results 
that one would obtain by solving numerically the Hamiltonian eigenvalue equation for 
the ground and first-excited states (more on this later). The real beauty of it all is that 
analytic expressions can be obtained for all these values. But some of them, however, fail 
to approximate the exact quantum results. Notably, the fluctuations of the quadrature 
operator Aq and of the atomic transition operator AJ X , obtained from the CS, show 
an erroneous behaviour. The reason can be traced down to the fact that, while the 
Hamiltonian ^ shows a parity symmetry given by 

[e i7rA , H] = (16) 

the coherent states mix states with opposite symmetry. 



3. Symmetry- Adapted States 

Considering the unitary transformations U (x) = exp (i x A.) , and using the relations 

U(x)aW( X ) = e-'xa, U( X ) J + U\ X ) = ^ x J+, we find 

U (x) H U\x) = o f a + oj a J z + -JL (a) J_ + a J + ) 

+ (e~ 2ix a) J + + e 2ix a J_) (17) 



'N 

so we have a symmetry transformation for x — 0, 7r. I.e., H is invariant under the group 
C 2 = {/,e^ A }. 

This symmetry allows for the classification of the Hamiltonian eigenstates in terms 
of the parity of the eigenvalues A = v + m + j of A, whereas states with opposite 
symmetry are strongly mixed by the variational CS states considered in the last section. 
We may, however, obtain symmetry- adapted states (SAS) [2H [22] by projecting with 
P± = \ (l±e i7rA ), i.e., 

k C)± =A4 P±\a, 

= A4 [ I a) ® |C) ± I - a) ® I - C) ] (18) 

with J\[± the normalisation factors for the even (+) and odd (— ) states. The surface 
energy takes a more complicated form 

2 1 



2^ H 1 \ l±e±(p 2 +i 2 )(cos9) 
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Figure 4. (Colour online.) Left: as function of the interaction strength 7, 

for the ground (blue, discontinuous upper curve) and first-excited (red, continuous 
lower curve) SAS states. Here, N — 20, to a — 1, and (j) c = 0. Right: q c /\/N as a 
function of 9 C for the CS and odd SAS states (continuous, green) and even SAS state 
(discontinuous, blue). 




... +1 tan 2 9 cos 9 

cos 9) ± 



l±e ± (P a +9 3 )(cos0) : F^J 
±p tan 9 sin + q e p2+q2 sin cos (cos #) _7V 



(19) 



eP 2 +9 2 (cos^)- w ± 1 

and the critical points which minimise it can be obtained, this time numerically, following 
the same method as for the coherent states (CS). In this case p c = 0; C = 0, 7r; and 
q c and 9 C now depend strongly on N. [Note: one could still obtain analytic expressions 
if, as usually done in many-body physics, one were to use the critical points for the 
CS energy surface, in order to obtain the trial state which approximates the two lowest 
symmetry-adapted energy states, and in which to evaluate the expectation values of 
the observables. Around the phase transition there would be observable differences, 
however, and one should calculate the critical points of (19) as we do here. For the first 
approach see 



Figure [4] (left) shows the behaviour of q c vs. 7 for the ground (even, blue, 
discontinuous upper curve) and first-excited (odd, red, continuous lower curve) states. 
The discontinuity in the phase space variable for the even state marks the phase 
transition [27]. This discontinuity arises from a competition for the global minimum 
in the even energy surface between two local minima, as 7 crosses the phase transition 
region, which can be seen in Fig. [5] [28]. In the case of the odd states, an inflection 
point marks the transition from one regime to the other. 

If we compare the minimum points q c associated to the energy surfaces of the CS 
and the even SAS, the discontinuity mentioned above translates into a forbidden range 
of values for 9 C (0.3 < 9 < 0.6) when q c is plotted against it in Fig. [4] (right). But it is 
interesting to note that, in the allowed ^-region, both fall along the same curve. The 
odd SAS state is also plotted, falls along the same curve, has no forbidden regions, and 
is indistinguishable from that of the coherent state. 
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7 = 0.545 7 = 0-552 & r = 0.560 




Figure 5. The symmetry adapted energy surface showing the minima for 7 = 
0.545,0.552,0.560 plotted as a function of q. Same values for N and uja as before. 



Finally, it must be mentioned that the SAS states have been successfully used to 
study phase transitions through the Renyi-Wehrl entropy [29] . 

4. Exact Solution 

In order to study exact quantum solutions we propose the correspondence q c — > 
=F ^2(0^0) , 6 C — > arccos ( ~ j J ^ ) 7 where =p corresponds to C = 0, tt. As the system is 
not integrable, we solve numerically the Hamiltonian eigenvalue equation for the ground 
and first-excited states, separating the even- and odd-parity cases. The quantum phase 
transition may be found through the fidelity between neighbouring states [27] [30] 

T = 1(^)1 Hi + 5 1 ))\ 2 , 

which at the transition acquires its minimum value. An associated, more precise method, 
is to use the fidelity susceptibility (the second derivative of the fidelity function given 
above) which is maximum at the phase transition region [31] . The dependence with 7 of 
the quadrature (or expectation value of the number operator) is shown in Figure [6] (left), 
for iV = 20 (bottom continuous gray dotted line) and N = 60 (top continuous black 
dotted line). These are compared with those for the SAS state at the same values of N 
(discontinuous dotted lines) and with the result for the CS (top continuous line). Note 
that the quantum and SAS solutions agree at (N ph ) = q c = 0, and approximate very 
well beyond the phase transition, while the CS fails to reproduce results well enough in 
the normal region. Recall that the SAS curves are discontinuous at the phase transition 
due to a competition for the global minimum in the even energy surface between two 
local minima, as explained in the last Section. After the discontinuity, the curves for 
any value of N become indistinguishable with the CS curve, as the figure shows. As 
N increases, 7 C — > 0.5 and the jump in the SAS curves becomes smaller; these curves 
approach more and more their quantum counterparts. Also, for greater values of 7 all 
curves tend to one and the same curve. 

One can also plot q c vs. 6 C (or \J (ala)/j vs. arccos((— J z )/j)) for the quantum 
state (Fig. [6] (right)) falling indistiguishably on the same universal curve together with 
the CS and SAS even and odd states. 
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Figure 6. (Colour online.) Left: dependence with 7 of the quadrature (expectation 
value of the number operator) for N = 20 (bottom continuous gray dotted line), and 
N = 60 (top continuous black dotted line), compared with those for the SAS state at 
the same values of N (discontinuous dotted lines) and with the result for the CS (top 
continuous cyan line). Right: Universal curve for quantum, SAS, and CS states. 



5. Critical Exponents 

A function f(r) is called homogeneous if, for all values of A, we have f(Xr) = g(X) f(r). 
The function g(X) is called the scaling function and must be of the form g(X) = X s ; s 
is called the critical exponent. The scaling hypothesis states that the singular part of 
many potentials are homogeneous functions near second-order phase transitions; this is 
true, for instance, for all thermodynamic potentials [32] . 

Critical exponents thus serve to describe the behaviour of important observables 
of a system near phase transitions and, although in most cases their properties are 
only supported by experimental data, they are believed to be universal with respect to 
physical systems. 

Our treatment for finite 2-level systems in a cavity, in the presence of a radiation 
field, allows us to study the critical value of the atom-field coupling parameter j c as a 
function of the number of atoms N, from which its critical exponent may be derived. 
Figure [7] shows this relationship for the exact, quantum ground (left) and first-excited 
(right) states. The points correspond to the numerical solution obtained from the 
diagonalisation of the Hamiltonian, while the continuous curves are the fits 

1 1 _2 

7c = g + 3 3 3 : (ground state) (20) 

1 2 _2 

7 C = — I — j^ 5 : (first-excited state) (21) 

2 5 

The value of j (= N/2) was varied from 5 to 900 for the ground state, and from 5 
to 450 for the first excited state. We see that the critical exponent is the same for both 
curves. A more demanding fit is the logarithmic function of both variables. When we 
take j G [100, 900] for the ground state we get the log-log plot shown in Fig. |8} Our 
best fit is 

ln[ 7c --] = -1.11596 ln[j] , (22) 



Mathematical Methods in Quantum Optics: the Dicke Model 



12 



Ground State 



First Excited State 



0.60 
0.58 
0.56 : 
0.54 : 
0.52 : 



200 



7, 
0.64. 

0.62 

0.60 

0.58 

0.56 

0.54 

0.52 



400 



600 



800 



200 



300 



400 



Figure 7. Behaviour of the critical value of the coupling parameter 7 C with the number 
of atoms N (= 2j), for the quantum ground (left) and first-excited (right) states. In 
both cases 7 C grows as N~ 2 / 3 . 
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Figure 8. Logarithmic behaviour of the critical value of the coupling parameter j c 
with the number of atoms N (= 2j), for the quantum ground state. The critical 
exponent is —2/3. 



and we note that e^ 1 ' 11596 = 0.32 76 ~ 1/3. 

While the SAS solutions do approximate very well the quantum solutions, except in 
a small vicinity of the phase transition, when we evaluate the critical coupling parameter 
coming from SAS as a function of the number of particles, we get a critical exponent of 
— 11/21, close but not exactly the appropriate one. In order to obtain this result, recall 
that the discontinuity in the phase space variables for the even state marks the critical 
value 7 C sc where the phase transition takes place. We varied j (= N/2) from 10 to 500 
for the ground state. Figure [9] (left) shows the results together with the fit (continuous 
curve) 

1 1 ii 

7c SC =2 + gf 5 . (23) 

The log-log plot for j E [100, 500] shown in Fig. [9] (right) confirms this value. Our best 
fit is 

Hw ~\}= ln[^] - ^ lnfol . (24) 
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Figure 9. Left: behaviour of the critical value of the coupling parameter 7 CSC with the 
number of atoms N (= 2j), for the SAS ground state. Right: logarithmic behaviour, 
showing a critical exponent of —11/21. 



with standard deviation a = 0.0059 and confidence interval for the slope (at a confidence 
level of 0.95) m E [-0.5337, -0.5154]. 

6. Discussion and Conclusions 

We have shown how various mathematical formalisms, specifically the catastrophe 
formalism and the use of group theory, may aid in the study of non-integrable systems 
in quantum optics, by comprehensively treating the problem of matter in the presence 
of a radiation field in a cavity, as modeled by the complete Dicke Hamiltonian. 

We approximated the phase transition of the Dicke model for a finite number N of 
atoms, via 3 different methods: i) numerical diagonalisation of the Hamiltonian and the 
use of the fidelity between neighbouring states; ii) variational (CS) test states that are 
a direct product of coherent Heisenberg-Weyl (ifW / (l))-states (for the electromagnetic 
field), and 5*(7(2)-states (for the atomic field); and iii) the use of projection operators on 
the coherent states in (ii) to obtain states which obey the parity symmetry in the total 
excitation number, present in the Hamiltonian. We saw that the symmetry- adapted 
states (SAS) in (iii) allow for the description of the dependence of observables in the 
number of particles N for the exact quantum ground and first-excited states. 

The three different approaches mentioned above lead to universal parametric curves 
for the expectation value of the first quadrature of the electromagnetic field (expectation 
value of the number operator), as a function of the atomic relative population. These 
are valid for all values of the matter-field coupling parameter, and valid for both the 
ground and first-excited states. 

Finally, using these mathematical tools, the critical value of the atom-field coupling 
parameter 7 C was found as a function of the number of atoms N, from which its critical 
exponent was derived. This is useful to describe the behaviour of the system near its 
phase transition. 
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